Purpose: To measure the growth of each strain in media with pH adjusted using lactic acid to 3.5, 4, 4.5, 5, 5.5, 6, and 7 (no lactic acid)

Protocol:

library(tidyverse)
library(xlsx)
library(broom)
library(rstatix)
library(ape)
library(ggbeeswarm)
library(formattable)
library(ggpubr)
library(cowplot)
library(corrplot)
library(kableExtra)

library(grid)
library(gridGraphics)

`%!in%` <- negate(`%in%`)

Short cuts and file paths

strains <- c("14018", "315-A", "C0011E4", "JCP7275", "C0084H9", "JCP8108", "JCP8017B", "C0040C2" , "UM35",  "JCP8066", "C0093B3", "AMD", "C0096A1", "C0102A2",  "UM224", "Gv5-1", "C0179B4", "C0056B5", "Gv101", "C0100B2", "C0101A1", "CMW7778B")

strains0 <- c("14018", "315-A", "C0011E4", "JCP7275", "C0084H9", "JCP8108", "JCP8017B", "C0040C2" , "UM35",  "JCP8066", "C0093B3", "AMD", "C0096A1", "C0102A2",  "UM224", "5-1", "C0179B4", "C0056B5", "101", "C0100B2", "C0101A1", "CMW7778B")

species <- c("Gardnerella vaginalis", "Gardnerella sp. 2", "Gardnerella sp. 3", "Gardnerella piotii", "Gardnerella leopoldii", "Gardnerella swidsinskii", "Gardnerella sp. 7", "Gardnerella sp. 8", "Gardnerella sp. 10", "Gardnerella sp. 11", "Gardnerella sp. 12")


figureOut <- "../../experiments_figures"
suppTableOut <- "../../experiments_figures/supplementary_posthoc_tables.xlsx"

Import data

midlog_raw_data <- read_csv("./data/pH_midlog_raw_data.csv")
pH_raw_data <- read_csv("./data/pH_raw_data.csv")

Plot settings

theme_set(theme_bw())

speciesColor <- scale_color_manual(values=c("#CC79A7", "#D55E00", "#56B4E9", "#0072B2", "#E69F00", "#009E73", 
                                                 "darkgreen", "gold3", "darkred", "darkblue",
                                                 "mediumpurple4"))

speciesFill <- scale_fill_manual(values=c("#CC79A7", "#D55E00", "#56B4E9", "#0072B2", "#E69F00", "#009E73", 
                                                 "darkgreen", "gold3", "darkred", "darkblue",
                                                 "mediumpurple4"))

Set up growth

Propagate twice and then incubate to mid log. 1:10 Dilution for mid Log: Protocol 1: 4 hour incubation Protocol 2: 7 hour incubation Protocol 3: 15 hour incubation ## Clean data

setupDF <- midlog_raw_data %>%
  group_by(Number, BioRep, Species, Strain, Batch, Step) %>%
  mutate(OD=OD-mean(.$OD[.$Number=="Blank"]),
                      OD=case_when(OD>0~OD,
                      OD<=0~0)) %>%
  filter(Number!="Blank") %>%
  mutate(Step=factor(Step, levels=c("Overnight 1", "Overnight 2", "Mid Log"))) %>%
  mutate(Strain=factor(Strain, levels=strains),
         Species=factor(Species, levels=species),
         Sample=paste(Strain, Batch, BioRep, sep="."),
         BatchNum=str_extract(Batch, "(?<=[A|B])[0-9]"),
         BatchRep=paste(BatchNum, BioRep, sep="."),
         StrainN=as.numeric(Strain),
         rectFill=case_when(StrainN %% 2 == 0 ~ "A",
                           StrainN %% 2 != 0 ~ "B")) %>%
  ungroup %>%
  with_groups(c("Sample", "Number", "BioRep", "Batch", "BatchNum", "BatchRep", "StrainN", "rectFill", "Species", "Strain", "Protocol",  "Step"), summarize, ODm=mean(OD), ODsd=sd(OD))

Assess mid log growth for removing poor growers with OD600 < 0.05

setupDF %>%
  filter(Step=="Mid Log") %>%
  ggplot(aes(x=Strain, y=ODm)) +
  geom_point(aes(shape=BioRep)) +
  geom_hline(aes(yintercept=0.05), linetype=2, color="gray") +
  theme(axis.text.x = element_text(angle=45, hjust=1)) +
  labs(x=NULL, y=bquote(OD[600]), shape="Biological Replicate", color="Experiment Batch")

failGrowthSamples <- setupDF %>%
  filter(Step=="Mid Log", ODm<0.05) %>%
  .$Sample
  
failGrowthSamples
##  [1] "C0100B2.1B1.A" "C0100B2.1B1.B" "C0100B2.1B2.A" "C0100B2.1B2.B"
##  [5] "C0102A2.1B1.A" "C0102A2.1B2.B" "C0179B4.1B1.A" "C0179B4.1B1.B"
##  [9] "C0179B4.1B2.A" "C0179B4.1B2.B" "Gv101.1B1.B"   "Gv5-1.1A2.A"  
## [13] "JCP7275.1B1.A" "JCP7275.1B1.B" "JCP7275.1B2.A" "JCP7275.1B2.B"
## [17] "JCP7275.1B3.A" "JCP7275.1B3.B"
overnightPlot <- setupDF %>%
  filter(Sample %!in% failGrowthSamples) %>%
  ggplot(aes(x=Strain, y=ODm)) +
  geom_rect(aes(xmin=StrainN-.5, xmax=StrainN+.5, ymin = 0, ymax = 1, fill=rectFill), alpha = 0.2, show.legend = FALSE) +
  geom_pointrange(aes(ymin=(ODm-ODsd), ymax=(ODm+ODsd), shape=BatchNum, color=Species), position = position_quasirandom(), size=0.3, alpha=0.7, show.legend = FALSE) +
  geom_point(aes(shape=BatchNum), alpha=0) +
  speciesColor +
  scale_fill_manual(values=c("white", "gray"), labels=c("A","B")) +
  ylim(0,1) +
  facet_grid(.~Step) +
  coord_flip() +
  scale_x_discrete(limits=rev) +
  theme(axis.text.x = element_text(size=8),
        axis.text.y = element_text(size=8),
        axis.title = element_text(size=11),
        panel.grid = element_blank(),
        legend.position = "none") +
  labs(x=NULL, y=bquote(OD[600]), shape="Experiment Batch")

## annotation portions of plot
#protocol bar
protocolBar <- setupDF %>%
  select(Strain, Species, Protocol) %>%
  unique %>%
  mutate(Protocol=as.character(Protocol)) %>%
  ggplot() +
    geom_bar(aes(x=Strain, y=2, fill=Protocol), stat="identity", width = 1.1) +
    scale_x_discrete(limits=rev) +
    scale_fill_manual(values=c("cornflowerblue", "burlywood1", "forestgreen")) +
    coord_flip() +
    theme_void() +
    theme(legend.position = "none")

speciesBar <- setupDF %>%
  select(Strain, Species, Protocol) %>%
  unique %>%
  ggplot() +
    geom_bar(aes(x=Strain, y=1, fill=Species), stat="identity", width=1.1) +
    scale_x_discrete(limits=rev) +
    speciesFill +
    coord_flip() +
    theme_void() +
    theme(legend.position = "none")

#Legends
protocol_legend <- ggpubr::get_legend(
    protocolBar + 
    guides(color = guide_legend(nrow = 1)) +
    theme(legend.position = "bottom"))
species_legend <- ggpubr::get_legend(
    speciesBar +
    guides(color=guide_legend(ncol = 1)) +
    theme(legend.position = "bottom"))
experiment_legend <- ggpubr::get_legend(
  overnightPlot +
  guides(shape=guide_legend(nrow = 1, override.aes = list(alpha=1, size=4))) +
  theme(legend.position = "bottom"))

a <- plot_grid(overnightPlot, speciesBar, protocolBar, nrow = 1, align = "h", axis = "bt", rel_widths = c(1, 0.02, 0.02))
b <- plot_grid(a, as_ggplot(species_legend), nrow=2, rel_heights = c(1, 0.3))
c <- plot_grid(as_ggplot(protocol_legend), as_ggplot(experiment_legend), nrow=1)
setUpGrowthFigure <- plot_grid(b, c, ncol=1, nrow=2, rel_heights = c(1, 0.1))
setUpGrowthFigure

Calculate growth in media with lactic acid-adjusted pH

Measured as OD600 in stationary phase ## Organize data

blankDF <- pH_raw_data %>%
    filter(Number=="Blank")

pH_results_raw <- pH_raw_data %>%
    filter(Number !="Blank")

Assess Blanks

blankDF %>%
  ggplot(aes(x=Plate, y=OD)) +
  geom_point(alpha=0.5) +
  facet_grid(pH~Batch, scales = "free_x")

May need to remove pH 6 data from protocol 3 strains.

Perform blank subtraction

meanBlankODs <- blankDF %>%
  with_groups(c(Plate, Batch, pH), summarize, blankOD=mean(OD))


DF0 <- pH_results_raw %>%
  left_join(meanBlankODs) %>%
  mutate(OD=OD-blankOD,
         Species=factor(Species, levels=species),
         Strain=factor(Strain, levels=strains),
         Sample=paste(Batch, Strain, BioRep, sep=".")) %>%
  filter(Sample %!in% failGrowthSamples)

Crude assessment of results

DF0 %>%
  ggplot(aes(x=pH, y=OD, color=pH, shape=BioRep)) +
  geom_point() +
  geom_hline(yintercept = 0, linetype=2, color="red") +
  facet_wrap(~Strain) +
  theme(axis.text.x = element_text(size=5)) +
  labs(x="pH", y=bquote(OD[600]), shape="Biological Replicate")

May be some issues with some pH = 7 samples.

Look at issues with low pH 7 results

DF0 %>%
  filter(Strain %in% c("14018", "315-A", "C0011E4", "AMD", "Gv5-1", "Gv101"), 
         pH==7) %>%
  mutate(PlateSet=str_extract(Batch, "(?<=[A|B])[1|2]"), 
         BatchRep=paste(Batch, BioRep, sep="."),
         BioRep=case_when(BatchRep=="1A1.A"~"A",
                          BatchRep=="1A1.B"~"B",
                          BatchRep=="1A2.A"~"C",
                          BatchRep=="1A2.B"~"D",
                          BatchRep=="1B1.A"~"E",
                          BatchRep=="1B1.B"~"F",
                          BatchRep=="2A1.A"~"A",
                          BatchRep=="2A1.B"~"B",
                          BatchRep=="2A2.A"~"C",
                          BatchRep=="2A2.B"~"D")) %>%
  ggplot(aes(x=Strain, y=OD, color=BioRep)) +
  geom_quasirandom(size=2, alpha=0.5) +
  labs(y=bquote(OD[600]), x=NULL)

All of Gv101 biological replicate A and one technical replicate of AMD

filterWells <- DF0 %>%
  filter(pH==7, 
         OD<0.005) %>%
  mutate(SampleWell=paste(Sample, well, sep=".")) %>%
  .$SampleWell
filterWells
##  [1] "1A1.Gv101.A.H7"   "1A1.Gv101.A.H8"   "1A1.Gv101.A.H9"   "1A2.Gv5-1.A.H1"  
##  [5] "1A2.Gv5-1.A.H2"   "1A2.Gv5-1.A.H3"   "1B1.C0102A2.A.H1" "1B1.C0102A2.A.H2"
##  [9] "1B1.Gv101.B.H4"   "2A2.AMD.A.H3"

Filter sample with no pH 7 values

filterSamples <- DF0 %>%
  mutate(SampleWell=paste(Sample, well, sep=".")) %>%
  filter(SampleWell %!in% filterWells) %>%
  count(Sample, pH) %>%
  spread(pH, n) %>%
  filter(is.na(`7`)) %>%
  .$Sample
filterSamples
## [1] "1A1.Gv101.A" "1A2.Gv5-1.A"

Filter pH 6 protocol 3

filterpHwells <- DF0 %>%
  mutate(SampleWell=paste(Sample, well, sep=".")) %>%
  filter(SampleWell %!in% filterWells,
         Protocol=="3",
         pH==6) %>%
  unique %>%
  .$SampleWell
filterpHwells
##  [1] "3A1.C0084H9.A.A6"  "3A1.C0084H9.A.B6"  "3A1.C0084H9.A.C6" 
##  [4] "3A1.C0084H9.B.D6"  "3A1.C0084H9.B.E6"  "3A1.C0084H9.B.F6" 
##  [7] "3A1.C0101A1.A.A12" "3A1.C0101A1.A.B12" "3A1.C0101A1.A.C12"
## [10] "3A1.C0101A1.B.D12" "3A1.C0101A1.B.E12" "3A1.C0101A1.B.F12"
## [13] "3A1.UM35.A.A6"     "3A1.UM35.A.B6"     "3A1.UM35.A.C6"    
## [16] "3A1.UM35.B.D6"     "3A1.UM35.B.E6"     "3A1.UM35.B.F6"    
## [19] "3A1.C0093B3.A.A12" "3A1.C0093B3.A.B12" "3A1.C0093B3.A.C12"
## [22] "3A1.C0093B3.B.D12" "3A1.C0093B3.B.E12" "3A1.C0093B3.B.F12"
## [25] "3A2.C0084H9.A.A6"  "3A2.C0084H9.A.B6"  "3A2.C0084H9.A.C6" 
## [28] "3A2.C0084H9.B.D6"  "3A2.C0084H9.B.E6"  "3A2.C0084H9.B.F6" 
## [31] "3A2.C0101A1.A.A12" "3A2.C0101A1.A.B12" "3A2.C0101A1.A.C12"
## [34] "3A2.C0101A1.B.D12" "3A2.C0101A1.B.E12" "3A2.C0101A1.B.F12"
## [37] "3A2.UM35.A.A6"     "3A2.UM35.A.B6"     "3A2.UM35.A.C6"    
## [40] "3A2.UM35.B.D6"     "3A2.UM35.B.E6"     "3A2.UM35.B.F6"    
## [43] "3A2.C0093B3.A.A12" "3A2.C0093B3.A.B12" "3A2.C0093B3.A.C12"
## [46] "3A2.C0093B3.B.D12" "3A2.C0093B3.B.E12" "3A2.C0093B3.B.F12"

Create final data frame

DF <- DF0 %>%
  mutate(SampleWell=paste(Sample, well, sep=".")) %>%
  filter(SampleWell %!in% filterWells,
         Sample %!in% filterSamples,
         SampleWell %!in% filterpHwells)

Plot

annotationDF <- DF %>%
  filter(pH==7) %>%
  with_groups(c(Number, Species, Strain, pH, Batch, Plate, BioRep, Protocol), summarize, mOD=mean(OD), sdOD=sd(OD)) %>%
  with_groups(c(Species, Strain), summarize, y=mean(mOD)) %>%
  with_groups(Species, mutate, y=1.2-.2*rank(Strain))

pHGrowthPlot <- DF %>%
  with_groups(c(Number, Species, Strain, pH, Batch, Plate, BioRep, Protocol), summarize, mOD=mean(OD), sdOD=sd(OD)) %>%
  mutate(pH=as.numeric(pH)) %>%
  ggplot() +
  stat_summary(aes(x=pH, y=mOD, group=Strain, linetype=Strain, color=Strain), fun="mean", geom = "line", show.legend = FALSE) +
  geom_pointrange(aes(x=pH, y=mOD, ymin=(mOD-sdOD), ymax=(mOD+sdOD), color=Strain, shape=Strain),size=0.3, show.legend = FALSE) +
  geom_text(data=annotationDF, aes(label=Strain, color=Strain, x=3.5, y=y), size=3, show.legend = FALSE, hjust=0) +
  facet_wrap(~Species) +
  scale_shape_manual(values=c(15, 16, 17, 15, 16, 17, 15, 16, 15, 16, 17, 15, 16, 15, 16, 17, 15, 15, 16, 15, 15, 15)) +
  scale_linetype_manual(values=c(1, 2, 3, 1, 2, 3, 1, 2, 1, 2, 3, 1, 2, 3, 1, 2, 1, 1, 2, 1, 1, 1)) +
  scale_color_manual(values=c("black", "gray47", "brown", "black", "gray47", "brown", "black", "gray47", "black", "gray47", "brown", "black", "gray47", "black", "gray47",  "brown", "black", "black", "gray47", "black", "black", "black")) +
  #theme(panel.grid = element_blank()) +
  labs(y=bquote(OD[600]), x="pH")
pHGrowthPlot

Minimum pH with non-zero growth per strain

T-test

ttest <- DF %>%
  with_groups(c(Number, Species, Strain, pH, Batch, Plate, BioRep, Protocol), summarize, mOD=mean(OD), sdOD=sd(OD)) %>%
  group_by(Species, Strain, pH) %>%
  t_test(mOD~1, mu=0.01, alternative = "greater") %>%
  add_significance("p")
minpHtable <- ttest %>%
  group_by(Strain, Species) %>%
  filter(p.signif!="ns") %>%
  filter(pH==min(pH)) %>%
  select(Species, Strain, pH) %>%
  arrange(pH) %>%
  kbl(caption="Table 1. Minimum pH with growth") %>%
  kable_classic(full_width=TRUE, html_font = "Arial")

minpHtable
Table 1. Minimum pH with growth
Species Strain pH
Gardnerella vaginalis 315-A 3.5
Gardnerella vaginalis C0011E4 4.5
Gardnerella sp. 2 C0084H9 4.5
Gardnerella piotii C0093B3 4.5
Gardnerella swidsinskii Gv5-1 4.5
Gardnerella sp. 10 C0100B2 4.5
Gardnerella sp. 11 C0101A1 4.5
Gardnerella vaginalis 14018 5.0
Gardnerella sp. 2 JCP7275 5.0
Gardnerella sp. 2 JCP8108 5.0
Gardnerella sp. 3 JCP8017B 5.0
Gardnerella piotii JCP8066 5.0
Gardnerella leopoldii AMD 5.0
Gardnerella leopoldii C0096A1 5.0
Gardnerella swidsinskii C0102A2 5.0
Gardnerella swidsinskii UM224 5.0
Gardnerella sp. 7 C0179B4 5.0
Gardnerella sp. 8 Gv101 5.0
Gardnerella sp. 3 C0040C2 5.5
Gardnerella piotii UM35 5.5
Gardnerella sp. 8 C0056B5 5.5
Gardnerella sp. 12 CMW7778B 5.5 
#save_kable(minpHtable, file = file.path(figureOut, paste(Sys.Date(), "Table1_minpHTable.png", sep="_")), zoom=6)

Moran’s I

Is minimum pH phylogenetically correlated by Moran’s I?

Load tree

We choose the weights as \(_{wij} = 1/d_{ij}\) , where the d’s is the distances measured on the tree:

We can now perform the analysis with Moran’s I:

Compare growth at pHs

Isolate Anova

One-way Anovas

DF %>%
  with_groups(c(Number, Species, Strain, pH, Batch, Plate, BioRep, Protocol), summarize, mOD=mean(OD), sdOD=sd(OD)) %>%
  select(Species, Strain, Protocol, pH, mOD, sdOD) %>%
  group_by(pH) %>%
  anova_test(mOD~Strain) %>%
  formattable()
pH Effect DFn DFd F p p<.05 ges
3.5 Strain 21 70 4.098 4.37e-06
0.551
4.0 Strain 21 70 1.999 1.70e-02
0.375
4.5 Strain 21 70 13.435 5.49e-17
0.801
5.0 Strain 21 70 23.869 5.59e-24
0.877
5.5 Strain 21 70 5.071 1.32e-07
0.603
6.0 Strain 17 58 2.217 1.30e-02
0.394
7.0 Strain 21 70 4.795 3.47e-07
0.590

Tukey’s HSD

odPostHocs <- DF %>%
  with_groups(c(Number, Species, Strain, pH, Batch, Plate, BioRep, Protocol), summarize, mOD=mean(OD), sdOD=sd(OD)) %>%
  select(Species, Strain, Protocol, pH, mOD, sdOD) %>%
  group_by(pH) %>%
  tukey_hsd(mOD~Strain, p.adjust.method = "BH")

# odPostHocs %>%
#  arrange(pH, p.adj) %>%
#  write.xlsx(file=suppTableOut, sheetName="TablesS7_pH_posthocs", col.names=TRUE, append=TRUE)
pHODpvals <- odPostHocs %>%
  split(.$pH) %>%
  map(~select(.x, group1, group2, p.adj)) %>%
  map(~dplyr::rename(.x, var1=group2, var2=group1, cor=p.adj)) %>%
  map(cor_spread) %>%
  map(~column_to_rownames(.x, var="rowname")) %>%
  map(as.matrix)

pHODcorplotInput <-  odPostHocs %>%
  split(.$pH) %>%
  map(~select(.x, group1, group2, estimate)) %>%
  map(~dplyr::rename(.x, var1=group2, var2=group1, cor=estimate)) %>%
  map(cor_spread) %>%
  map(~column_to_rownames(.x, var="rowname")) %>%
  map(as.matrix)

makeCorrPlots <- function(inputMatrix, pvalMatrix) {
  #dev.new()
  #par(xpd = NA, # switch off clipping, necessary to always see axis labels
  #bg = "transparent", # switch off background to avoid obscuring adjacent plots
  #oma = c(2, 2, 0, 0)) # move plot to the right and up
  corrplot(inputMatrix, is.corr=FALSE, type = "lower", tl.col="black", p.mat = pvalMatrix, insig = "label_sig", sig.level = c(0.001, 0.01, 0.05), pch.cex = 0.9)
  grid.echo()
  plot <- grid.grab()
  #dev.off()
  return(plot)
}

pHODcorrPlots <- map2(pHODcorplotInput, pHODpvals, makeCorrPlots)

Moran’s I

Proportion of pH 7 growth

ratioDF <- DF %>%
  with_groups(c(Number, Species, Strain, pH, Batch, Plate, BioRep, Protocol), summarize, mOD=mean(OD)) %>%
  spread(pH, mOD) %>%
  mutate_at(vars("3.5", "4", "4.5", "5", "5.5", "6", "7"), ~.x/`7`) %>%
  gather("pH", "ratio", c("3.5", "4", "4.5", "5", "5.5", "6", "7"))

Plot

annotationDFratio <- DF %>%
  filter(pH==7) %>%
  with_groups(c(Number, Species, Strain, pH, Batch, Plate, BioRep, Protocol), summarize, mOD=mean(OD), sdOD=sd(OD)) %>%
  with_groups(c(Species, Strain), summarize, y=mean(mOD)) %>%
  with_groups(Species, mutate, y=2.5-.5*rank(Strain))

pHratioPlot <- ratioDF %>%
  mutate(pH=as.numeric(pH)) %>%
  ggplot() +
  stat_summary(aes(x=pH, y=ratio, group=Strain, linetype=Strain, color=Strain), fun="mean", geom = "line", show.legend = FALSE) +
  geom_point(aes(x=pH, y=ratio, color=Strain, shape=Strain), show.legend = FALSE) +
  geom_text(data=annotationDFratio, aes(label=Strain, color=Strain, x=3.5, y=y), size=3, show.legend = FALSE, hjust=0) +
  facet_wrap(~Species) +
  scale_shape_manual(values=c(15, 16, 17, 15, 16, 17, 15, 16, 15, 16, 17, 15, 16, 15, 16, 17, 15, 15, 16, 15, 15, 15)) +
  scale_linetype_manual(values=c(1, 2, 3, 1, 2, 3, 1, 2, 1, 2, 3, 1, 2, 3, 1, 2, 1, 1, 2, 1, 1, 1)) +
  scale_color_manual(values=c("black", "gray47", "brown", "black", "gray47", "brown", "black", "gray47", "black", "gray47", "brown", "black", "gray47", "black", "gray47",  "brown", "black", "black", "gray47", "black", "black", "black")) +
  labs(y="Proportion pH 7 growth", x="pH")
pHratioPlot

Compare

Average proportion pH 7 growth at each pH

ratioDF %>%
  filter(pH!=7) %>%
  with_groups(pH, summarize, perc=mean(ratio, na.rm=TRUE)*100)
## # A tibble: 6 × 2
##   pH       perc
##   <chr>   <dbl>
## 1 3.5    1.79  
## 2 4      0.0708
## 3 4.5    4.40  
## 4 5     17.8   
## 5 5.5   37.4   
## 6 6     58.4

Isolate Anova

One-way Anovas

ratioDF %>%
  filter(pH!="7") %>%
  select(Species, Strain, Protocol, pH,ratio) %>%
  group_by(pH) %>%
  anova_test(ratio~Strain) %>%
  formattable()
pH Effect DFn DFd F p p<.05 ges
3.5 Strain 21 70 0.293 9.99e-01 0.081
4 Strain 21 70 1.129 3.41e-01 0.253
4.5 Strain 21 70 0.608 8.99e-01 0.154
5 Strain 21 70 2.005 1.60e-02
0.376
5.5 Strain 21 70 2.377 4.00e-03
0.416
6 Strain 17 58 5.565 3.87e-07
0.620

Tukey’s HSD

ratioPostHocs <- ratioDF %>%
  filter(pH!="7") %>%
  select(Species, Strain, Protocol, pH,ratio) %>%
  group_by(pH) %>%
  tukey_hsd(ratio~Strain, p.adjust.method = "BH")

# ratioPostHocs %>%
#  arrange(pH, p.adj) %>%
#  write.xlsx(file=suppTableOut, sheetName="TablesS8_pH_ratio_posthocs", col.names=TRUE, append=TRUE)
pHratiopvals <- ratioPostHocs %>%
  filter(pH!="7") %>%
  split(.$pH) %>%
  map(~select(.x, group1, group2, p.adj)) %>%
  map(~dplyr::rename(.x, var1=group2, var2=group1, cor=p.adj)) %>%
  map(cor_spread) %>%
  map(~column_to_rownames(.x, var="rowname")) %>%
  map(as.matrix)

pHratiocorplotInput <-  ratioPostHocs %>%
  filter(pH!="7") %>%
  split(.$pH) %>%
  map(~select(.x, group1, group2, estimate)) %>%
  map(~dplyr::rename(.x, var1=group2, var2=group1, cor=estimate)) %>%
  map(cor_spread) %>%
  map(~column_to_rownames(.x, var="rowname")) %>%
  map(as.matrix)

pHratiocorrPlots <- map2(pHratiocorplotInput, pHratiopvals, makeCorrPlots)

Figure 5: Growth at varying pH’s

plot_grid(pHGrowthPlot+theme(axis.title.y = element_text(size=18), strip.text.x = element_text(size=7)), pHratioPlot+theme(axis.title.y = element_text(size=18), strip.text.x = element_text(size=7)), labels = c("A", "B"), label_size = 15)

#ggsave(filename = file.path(figureOut, paste(Sys.Date(), "Figure5_pHGrowth.png", sep = "_")))

Figure S4: Post Hoc Test Results

ODPostHocPlots <- plot_grid(pHODcorrPlots$`3.5`, pHODcorrPlots$`4`, pHODcorrPlots$`4.5`, pHODcorrPlots$`5`, pHODcorrPlots$`5.5`, pHODcorrPlots$`6`, pHODcorrPlots$`7`, ncol=1,  labels = c("A", "C", "E", "G", "I", "K", "M"), label_size = 20)

ratioPostHocPlots <- plot_grid(pHratiocorrPlots$`3.5`, pHratiocorrPlots$`4`,pHratiocorrPlots$`4.5`, pHratiocorrPlots$`5`, pHratiocorrPlots$`5.5`, pHratiocorrPlots$`6`, ggplot+geom_blank(),  ncol=1, labels = c("B", "D", "F", "H",  "J",  "L", " "), label_size = 20)

plot_grid(ODPostHocPlots, ratioPostHocPlots, ncol=2, align = "hv") 

#ggsave(file.path(figureOut, paste(Sys.Date(), "FigureS4_pH_posthocs.png", sep = "_")))
sessionInfo()
## R version 4.4.0 (2024-04-24)
## Platform: x86_64-apple-darwin20
## Running under: macOS Ventura 13.6.6
## 
## Matrix products: default
## BLAS:   /Library/Frameworks/R.framework/Versions/4.4-x86_64/Resources/lib/libRblas.0.dylib 
## LAPACK: /Library/Frameworks/R.framework/Versions/4.4-x86_64/Resources/lib/libRlapack.dylib;  LAPACK version 3.12.0
## 
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
## 
## time zone: America/New_York
## tzcode source: internal
## 
## attached base packages:
## [1] grid      stats     graphics  grDevices utils     datasets  methods  
## [8] base     
## 
## other attached packages:
##  [1] gridGraphics_0.5-1 kableExtra_1.4.0   corrplot_0.92      cowplot_1.1.3     
##  [5] ggpubr_0.6.0       formattable_0.2.1  ggbeeswarm_0.7.2   ape_5.8           
##  [9] rstatix_0.7.2      broom_1.0.6        xlsx_0.6.5         lubridate_1.9.3   
## [13] forcats_1.0.0      stringr_1.5.1      dplyr_1.1.4        purrr_1.0.2       
## [17] readr_2.1.5        tidyr_1.3.1        tibble_3.2.1       ggplot2_3.5.1     
## [21] tidyverse_2.0.0   
## 
## loaded via a namespace (and not attached):
##  [1] gtable_0.3.5      beeswarm_0.4.0    xfun_0.46         bslib_0.7.0      
##  [5] htmlwidgets_1.6.4 rJava_1.0-11      lattice_0.22-6    tzdb_0.4.0       
##  [9] vctrs_0.6.5       tools_4.4.0       generics_0.1.3    parallel_4.4.0   
## [13] fansi_1.0.6       highr_0.11        pkgconfig_2.0.3   lifecycle_1.0.4  
## [17] farver_2.1.2      compiler_4.4.0    munsell_0.5.1     carData_3.0-5    
## [21] vipor_0.4.7       htmltools_0.5.8.1 sass_0.4.9        yaml_2.3.9       
## [25] crayon_1.5.3      pillar_1.9.0      car_3.1-2         jquerylib_0.1.4  
## [29] cachem_1.1.0      abind_1.4-5       nlme_3.1-165      tidyselect_1.2.1 
## [33] digest_0.6.36     stringi_1.8.4     labeling_0.4.3    fastmap_1.2.0    
## [37] colorspace_2.1-0  cli_3.6.3         magrittr_2.0.3    xlsxjars_0.6.1   
## [41] utf8_1.2.4        withr_3.0.0       scales_1.3.0      backports_1.5.0  
## [45] bit64_4.0.5       timechange_0.3.0  rmarkdown_2.27    bit_4.0.5        
## [49] ggsignif_0.6.4    hms_1.1.3         evaluate_0.24.0   knitr_1.48       
## [53] viridisLite_0.4.2 rlang_1.1.4       Rcpp_1.0.13       glue_1.7.0       
## [57] xml2_1.3.6        vroom_1.6.5       svglite_2.1.3     rstudioapi_0.16.0
## [61] jsonlite_1.8.8    R6_2.5.1          systemfonts_1.1.0